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was  performed  under  In-house  Project  No.  241803,  "Dynamic  Behavior  of 
Engine  Materials."  The  work  was  conducted  by  Mr.  Jeffry  D.  Sharp, 
collocated  from  ASD/ENFSM,  Wrlght-Patterson  AFB,  Ohio. 

The  author  would  like  to  offer  thanks  to  his  thesis  advisor. 

Dr.  George  Sutherland,  Ohio  State  University,  for  his  suggestions  and 
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making  available  the  experimental  facilities  and  providing  technical  advice, 
and  to  Ms.  Barbara  Lear,  ASD/ENF  for  her  secretarial  assistance. 

The  research  was  conducted  during  the  period  October  1977  to 
July  1979.  This  report  was  submitted  for  publication  in  September  1979. 
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SUMMARY 

With  thrust- to-weight  ratios  in  aircraft  gas  turbine  engines 
increasing,  all  components  are  required  to  operate  at  higher  steady- 
stress  levels.  This  situation,  coupled  with  cyclic  life  limits  based 
on  crack  propagation,  demands  a  better  understanding  of  the  operating 
environment  and  various  loading  conditions  to  which  a  component  may  be 
subjected  throughout  Its  useful  life.  Vibratory  stresses  induced  in 
rotating  airfoils  by  foreign  object  impacts  are  a  real  part  of  jet 
engine  operating  environments  and  require  careful  consideration  when 
designing  high  performance  blades.  However,  before  one  can  analyze 
a  component  for  such  conditions,  the  loading  mechanism  must  be  under¬ 
stood.  This  has  been  investigated  for  members  of  similar  stiffnesses 
impacting  each  other  (References  1  through  4),  but  little  work  has  been 
done  to  define  and  analyze  the  effect  of  a  soft  body  impacting  a  member 
with  the  approximate  geometry  and  stiffness  of  a  fan  blade.  The  first 
step  toward  accomplishing  this  understanding  Is  to  fabricate  several 
cantilever  beam  test  specimens,  impact  them  with  soft  projectiles,  and 
record  the  strain  response  at  various  locations  as  a  function  of  time. 

From  this  data  base,  a  simple  beam  theory  model  can  be  evaluated  for 
several  different  cases,  each  treating  the  impact  as  a  slightly  dif¬ 
ferent  phenomenon,  but  based  on  the  common  assumption  that  all  the 
projectile  momentum  Is  transferred  to  the  beam.  The  results  indicate 
that  a  soft-body  Impact  on  a  cantilever  beam  can  be  modeled  as  forced 
vibration  for  the  duration  of  the  Impact,  then  free  vibration  for  all 
time  thereafter.  The  forcing  function  created  by  the  projectile  is 
best  modeled  as  a  step-function  distributed  over  an  area.  Damping  in 
the  system  has  a  significant  effect  on  strain  In  the  beam,  but  is  dif¬ 
ficult  to  accurately  predict,  and  thus  should  not  be  considered  in 
maximum  strain  calculations.  However,  an  accurate  strain/time  history 
can  be  predicted  if  modal  damping  values  for  the  first  four  resonant 
modes  are  known.  Impact  testing  of  very  large  structures  can  be  performed 
on  scaled  down  models  with  no  change  in  results. 
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SECTION  I 
INTRODUCTION 

The  problem  of  predicting  the  response  of  a  blade-like  structure  to 
an  Impact  by  a  soft  body  Is  very  complex  and  Involves  many  structural 
interactions.  To  begin  analyzing  all  the  load  transfer  mechanisms 
present  In  such  a  problem  Is  a  monumental  task.  One  must  initially  take 
a  more  basic  approach.  The  objective  of  this  study  Is  to  determine  how 
a  simple  cantilever  beam  responds  to  an  Impact  on  its  centerline  by  a 
soft  object.  If  the  overall  bending  strain  at  various  spanwise  locations 
in  a  beam  can  be  recorded  throughout  an  Impact,  the  gross  loading 
mechanism  can  be  studied  and  an  analytical  model  can  be  developed  to 
predict  the  response. 

To  accomplish  this  objective,  four  geometrically  similar  beams  were 
fabricated  and  instrumented  with  strain  gages  at  various  spanwise  loca¬ 
tions,  mounted  in  a  cantilever  configuration,  and  Impacted  with  geomet¬ 
rically  similar  soft  projectiles.  Data  from  the  strain  gages  was 
digitally  recorded  during  impact,  stored,  then  analyzed  to  determine 
maximum  strain  at  each  gage  location  and  the  modal  content  of  each 
response.  Any  damping  effects  present  were  noted.  This  provided  the 
dynamic  response  data  necessary  to  evaluate  several  proposed  analytical 
models. 

Next,  analytical  models  based  on  the  Euler-Bernoulli  Beam  Theory  and 
employing  linear  modal  analysis  were  formulated.  Two  basic  approaches 
were  used  to  simulate  the  Impact  loading.  The  first  was  to  assume  a 
purely  impulsive  loading  In  which  the  projectile  Imparted  an  instan¬ 
taneous  initial  velocity  to  the  beam.  The  second  treated  the  impact  as 
a  forced  vibration  problem,  which  modeled  the  impact  as  a.)  a  step 
function  and  b.)  a  half  sine  wave  In  time.  Each  model  was  defined  and 
formulated  on  a  computer  and  the  resulting  time  responses  were  compared 
with  the  test  data  to  determine  Its  accuracy.  The  following  text 
describes  in  detail  the  experimental  work  performed,  the  data  analysis 
techniques  used,  the  different  model  derivations  and  formulations,  and 
finally,  the  results  obtained  from  each  model  and  a  comparison  with 
experimental  results. 
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SECTION  II 

EXPERIMENTAL  PROCEDURE  AND  RESULTS 

Four  cantilever  beam  specimens  were  fabricated  from  7075  T-6 
aluminum,  approximately  91  R  or  190  B  hardness,  the  smallest  beam  being 
5.6  inches  long,  1.2  Inches  wide,  and  0.15  inches  thick.  The  remaining 
beams  were  scaled  from  these  dimensions  by  factors  of  1.5,  2.0,  and  2.5. 

To  ensure  cantilever  boundary  conditions  and  minimize  damping  effects 
from  mounting  fixtures,  the  root  section  of  each  specimen  was  left 
significantly  thicker  than  the  rest  of  the  beam,  as  shown  in  Figure  1. 

For  simplicity,  the  Impact  site  of  each  beam  was  chosen  to  be  at 
75%  span.  This  provided  a  reasonable  portion  on  either  side  of  the 
desired  impact  location  to  ensure  contact  by  the  whole  projectile. 

Strain  gages  were  mounted  on  each  beam  at  the  root  and  75%  span  on  the 
side  opposite  that  to  be  Impacted.  On  three  of  the  four  beams,  gages 
were  also  placed  at  25%  and  50%  span  to  record  additional  data  (Figure  2). 
These  four  locations  were  chosen  to  help  identify  the  location  of 
maximum  strain  and  to  provide  a  sufficient  amount  of  strain/time  history 
information  to  facilitate  verification  of  analytical  models. 

The  projectiles  used  to  impact  each  beam  were  spherical  bodies  of 
micro-balloon  gelatin,  a  porous  gelatin  used  to  simulate  bird  impacts 
(Reference  7).  Projectile  sizes  and  velocities  were  selected  to  meet 
several  criteria.  The  base  line  projectile  diameter  was  chosen  to  be 
0.5  inches  so  it  would  be  small  compared  to  the  smallest  beam  length 
and  width.  The  geometric  ratios  used  to  size  the  beams  were  also  used 
to  determine  the  remaining  projectile  diameters.  This  calculation 
yielded  projectile  diameters  of  0.5,  0.75,  1.0,  and  1.25  inches. 

However,  a  0.75  Inch  diameter  mold  was  not  available,  so  the  projectile 
sizes  actually  used  for  the  Impact  experiments  were  0.5,  0.70,  1.0,  and 
1.25  inches  in  diameter.  A  calculation  was  made  for  approximating 
projectile  velocities  based  on  a  one  degree-of-freedom  lumped  mass 
beam  model  (Reference  6).  To  assure  elastic  response,  velocity  calcula¬ 
tions  were  based  on  0.25%  strain  at  the  beam  root,  which  predicted  a 
projectile  velocity  of  about  570  fps.  These  calculations  are  detailed 
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in  Appendix  A.  To  account  for  the  gross  model  used,  the  velocity 
prediction  was  decreased  to  400  fps.  An  attempt  was  made  to  maintain 
this  projectile  velocity  for  each  beam  specimen  so  the  concept  of 
scaling  impact  specimens  could  be  Investigated. 

Two  different  devices  were  used  to  record  transient  strain  data. 

A  Hewlett-Packard  5451 B  Fast  Fourier  Analyzer  (FFT)  capable  of 
digitizing  and  storing  four  channels  of  data  at  20  KHZ  was  employed  to 
provide  a  record  of  transient  response  data  which  could  be  analyzed  at 
a  later  time.  A  Zonics  AE  102-2  transient  recorder  capable  of  digitizing 
eight  channels  of  data  from  2  KHZ  to  200  KHZ  was  utilized  to  Investigate 
short  time  and  long  time  responses.  Strain  data  recorded  by  the  FFT  unit 
was  stored  on  a  magnetic  disk,  while  the  transient  recorder  data  was 
stored  temporarily  In  the  recorder  and  then  input  to  an  oscilloscope  and 
photographed. 

For  testing,  each  beam  was  mounted  In  a  steel  fixture  as  shown  in 
Figure  3.  The  entire  assembly  was  then  placed  inside  a  steel  enclosure 
on  the  impact  range ,  positioned  with  a  laser  aiming  device  and  then 
secured  to  the  enclosure. 

Projectiles  were  first  weighed  and  the  weight  recorded,  and  then 
launched  from  a  smooth-bore  tube,  propelled  by  compressed  air,  compressed 
helium,  or  burning  gunpowder.  Each  projectile  was  carried  down  the 
tube  in  a  sabot,  a  plastic  bore  fitting  carrier,  which  protected  it 
from  the  launch  tube  walls.  Several  Inches  In  front  of  the  target  a 
constriction  in  the  barrel  stopped  the  sabot,  allowing  the  projectile 
to  continue.  Before  striking  the  beam,  the  projectile  tripped  a  pair  of 
laser  light  sources  connected  to  a  time  Interval  counter  to  measure  Its 
velocity. 

Several  shots  were  made  at  each  beam  to  "zero  in"  on  the  amount  of 
pressure/powder  necessary  to  attain  the  desired  velocity.  Impacts  on 
the  first  beam  Indicated  that  a  velocity  between  300  and  350  fps  was 
adequate  to  produce  the  desired  strain  at  the  specimen  root.  Thus, 
those  values  became  bounds  on  the  desired  velocity  for  all  remaining  tests. 
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Figures  4A  through  4D  are  strain/time  responses  from  each  beam  root 
location.  Variations  in  peak  strain  values  due  to  differing  projectile 
velocities  and  densities  can  be  seen.  The  data  indicates,  however,  that 
linear  scaling  of  all  beam  and  projectile  dimensions  will  produce  equal 
strains.  Appendix  B  is  a  dimensional  analysis  that  supports  this 
observation.  Responses  from  the  impact  site  (75%  span)  on  each  beam 
are  shown  in  Figures  5A  through  5C.  The  same  observations  made  of  the 
root  strain  data  can  also  be  made  of  these  traces.  In  addition,  very 
little  damping  effect  can  be  observed  in  either  set  of  data,  except  at 
high  frequency.  Figures  6A  through  6D  are  the  strain  responses  from  each 
gage  location  on  the  14.0  Inch  beam.  These  traces  indicate  that  the 
maximum  strain  experienced  during  impact  occurs  at  the  beam  root. 

Figures  7A  through  7D  show  the  response  of  the  14.0  inch  beam  over  0.5 
milliseconds  as  being  slow  and  containing  no  high  frequency  components. 
Figures  8A  through  80  are  the  responses  of  the  same  locations  over  500 
milliseconds,  demonstrating  damping  in  the  system  affecting  only  long¬ 
term  response.  All  the  strain  traces  Indicate  the  presence  of  several 
frequency  components.  However,  It  should  be  noted  that  the  sampling  rate 
used  to  record  this  data  was  too  slow  to  pick  up  frequencies  above  the 
first  or  second  mode.  Figures  9A  through  9D  are  Digital  Fourier  Trans¬ 
forms  of  root  strain  responses  from  each  beam,  showing  the  presence  of 
the  first  four  beam  bending  modes.  Similar  plots  for  the  75%  span 
location.  Figures  10A  through  10C  show  the  same  result  with  a  small 
contribution  from  the  fifth  mode.  Table  1  summarizes  the  maximum  strain 
values  from  all  experimental  records. 

An  overview  of  the  test  results  leads  to  several  conclusions.  The 
maximum  stress  over  the  entire  impact  event  occurs  at  the  beam  root. 

For  the  same  projectile  velocity  and  density,  linearly  scaled  beams  and 
projectiles  will  produce  the  same  impact  strain  values.  The  strain 
response  at  any  location  in  the  beam  has  significant  contribution  from 
only  the  first  five  resonant  bending  modes.  Damping  In  the  system  has 
little  effect  on  the  peak  strain  value  but  does  become  apparent  in  long 
time  response.  These  results  provided  the  basis  for  evaluating  several 
analytical  models. 
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Figure  5.  (Contd) 


♦No  Data  Available  At  This  Location  From  11.2"  Beam. 
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Experimental  Strain  in  14.0  Inch  Beam  (500  ms.) 
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Figure  9.  Fourier  Transforms  of  Root  Strain  Data 
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SECTION  III 
THEORETICAL  SOLUTION 

Two  basic  approaches  were  taken  to  modeling  the  soft-body  impact 
problem.  The  first  entailed  treating  the  projectile  as  imparting  an 
initial  velocity  to  the  beam  at  the  point  of  impact.  This  approach  was 
suggested  in  Reference  9  and  applied  to  a  hard-body  impact  problem  in 
Reference  8.  The  second  approach  was  to  treat  the  beam  response  to 
impact  as  a  forced  vibration,  as  in  Reference  4.  The  forcing  function 
was  modeled  as  both  a  square  wave  and  a  half-sine  wave.  In  both 
formulations,  the  common  assumption  made  was  that  all  the  projectile 
momentum  is  transferred  to  the  beam,  i.e.,  linear  momentum  is  conserved 
and  the  projectile  assumes  a  velocity  of  zero  after  Impact.  This 
assumption,  which  is  supported  in  Reference  7  and  by  experimental 
observation,  implies  that  the  coefficient  of  restitution  is  close  to 
zero,  but  unknown.  Both  formulations  are  based  on  linear  modal  theory 
and  employ  the  Euler-Bernoul 1 i  beam  relationships.  Timoshenko  beam 
effects  (i.e.,  shear  and  rotary  inertia)  were  investigated  and  found  to 
be  negligible  (Appendix  C). 

The  solution  for  free  vibration  of  a  beam  developed  in  Reference  5 
was  used  to  determine  the  resonant  frequencies  and  mode  shapes  of  a 
cantilever.  The  general  beam  mode  shape  (as  a  function  of  space 
coordinate  only)  is  given  by: 

y(x)  =  A  sinh(Bx)  +  B  cosh(Bx)  +  C  sin(3x)  +  D  cos(Bx)  (1  ) 

where  8  is  dependent  upon  the  applied  boundary  conditions.  For  a 
cantilever: 

y  *  0 

at  x  ■  0 

dy/dx  =  0 

M  =  0  or  d^y/dx^  =  0 

at  x  *  t 

V  -  0  or  d3y/dx3  -  0 
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Substitution  of  these  boundary  conditions  Into  Equation  1  yields: 

cosh(B£)cos(B£)  +  1  »  0  (2) 

which  is  the  solution  for  gb  for  the  various  resonant  modes.  The  mode 
shape  expression  can  also  be  rewritten,  as  a  result  of  applying  the 
boundary  conditions,  as: 

y(x)  =  A  (cosh(Bx)  -  cos(Bx)  -  K  [sinh(Bx)  -  sin(Bx)J)  (3) 

where  k  =  sinh(g£)  -  s1n(g£)/cosh(g«,)  +  cos(BJl).  For  simplicity. 

Equation  3  will  be  written  as: 

y(x)  «  A  (j>(x) 

defining 

<t» (x)  =  cosh(8x)  -  cos(8x)  -  K  [sinh(Bx)  -  sin(Bx)]  (3a) 

Employing  linear  modal  theory,  the  total  response  of  the  beam,  Including 
time  variation,  can  be  formulated  as: 

y(x,  t)  -  N  |  x  AN*H(x)TN(t)  (4) 

with  <J>N(x)  corresponding  to  a  solution  of  Equation  3a  for  a  particular 
value  of  8,  determined  from  Equation  2.  TN( t )  is  a  time  varying  function 
of  the  resonant  frequency  for  each  mode,  w,.,  calculated  from 

-  82n  EI/pA  (5) 

Assuming  a  sinusoidal  time  response, 

Tjj(t)  =  C1  sinCBjjt)  +  C2  cosO^t)  (6) 

and  C2  are  determined  by  evaluating  the  initial  conditions  of  the 
problem.  At  this  point,  each  approach  to  analyzing  the  soft  body  impact 
must  be  treated  separately. 

1.  INITIAL  VELOCITY  SOLUTION 

The  first  approach  taken  was  to  assume  the  projectile  imparts  an 
initial  velocity  to  the  beam,  as  shown  pictorially  in  Figure  11 A .  The 
velocity  can  be  expressed  as: 

V  -  V6(x  -  L)  (7) 


2 


y(x,  0)  «  V6(x  -  L) 


Application  of  the  I.C.s  yields  C2  =  0  and 

V  (x  -  L)  -  N  I  1  V^N(x)  (8) 

Multiplying  both  sides  by  a  particular  mode  shape,  ^(x),  applying 
orthogonality  conditions  and  integrating  over  the  length  of  the  beam 
results  in: 

'*V6(x  -  L)$M(x)dx  -  V^J^xldx  (9) 

The  left  hand  side  of  Equation  9  reduces  to  V4>M ( L ) »  while  the  integral 
on  the  right  hand  side,  when  evaluated,  equals  i  (Reference  10). 
Therefore, 

which  leads  to  the  solution  for  AN> 

■  A1(ta>1/u>N)  [d»N(L) /^^(L)  ]  (10) 

A1  must  now  be  solved  for  to  determine  the  complete  solution.  The  rela¬ 
tionship  remaining  is  conservation  of  momentum.  Applying  this  principle 
as  previously  described,  the  projectile  momentum  is  equated  to  the  net 
momentum  of  the  beam.  Thus, 

Vp  "  JOpbV(x*  t)dx 

Vp  "  ^  N  I  1  Vn  COS%t)  'J*ll(x)dx 

Evaluating  the  integral  and  substituting  Equation  10  for  A^,  we  have 

Vp  '  2pb*»  h  I  1  'Will'S. 
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where  is  defined  in  Equation  3a.  This  relationship  exists  only  at 
t  =  0,  therefore, 

Vp  =  2MBAl°'l/^l(L)  n  I  iVl>V®n 

Defining  C$  =  N  f  and  solving  for  A1 , 

Al  "  mpVp*i(l)/hbuics 

Thus,  one  can  solve  for  A^  by  substituting  into  Equation  10. 

*»-VtVl)Ws  l"> 

The  beam  response  to  impact  can  now  be  written  as: 

y(x,  c)  =  MpVp/f^Cg  N  |  1[*N(L)/wN]<(>N(x)cos(a)Nt)  (12) 

This  solution  can  be  further  modified  to  simulate  the  projectile 
impact  by  spreading  the  initial  velocity  out  from  a  point  to  an  area, 
shown  in  Figure  11 B.  Equation  7  is  rewritten  as: 

V  «  V(u[x  -  (L  -  Dp/2) ]  -  U[x  -  (L  +  Dp/2)]1  (13) 

which  spreads  the  impact  over  an  area  as  wide  as  the  projectile  and 
euqally  distributes  it  on  either  side  of  the  impact  span.  The  left  hand 
side  of  Equation  9  now  becomes: 

l*V(U[x  -  (L  -  Dp/2) )  -  U[x  -  (L  +  Dp/2)))«M(x)dx 

which  reduces  to: 


_  t  +  D  /2 
VJ  P  (k(x)dx 

■--V2 

Following  through  the  solution,  the  final  expression  for  the  beam 
response  becomes: 

y(x,  t)  =  MpVp/MBC's  N  |  1[41n/u1j)^(x)cos(^c)  (14) 
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where 

L  +  D  /2 

C  -J  P  $>(x)dx 
*  ^-Dp/2 


(15) 


and 

c's  - »  I  i2v\.  V  <’6> 

Both  solutions,  for  impact  at  a  point  and  over  an  area,  were  programmed 
on  the  computer  to  be  evaluated  against  test  data. 

2.  FORCED  VIBRATION 

The  second  approach  taken  was  to  treat  the  impact  as  a  forced 
vibration  problem  for  the  duration  of  projectile  contact,  shown  in 
Figure  12A.  The  solution  to  such  a  problem  is  in  two  parts,  the 
homogeneous  (free  vibration)  solution  and  the  particular  (forced  vibration) 
solution.  The  homogeneous  solution  is  that  given  in  Equation  4.  The 
particular  solution  is  developed  below. 

First,  the  forcing  function  was  assumed  to  be  a  constant  force  at 
a  point  acting  over  a  finite  time  period.  In  other  words, 

F(x,  c)  =  FS(x  -  L)[U(t)  -  U(t  -  T0)J  (17) 

The  time  interval,  Tq,  over  which  the  force  acts  is  determined  in 
Reference  7  to  be  the  time  necessary  for  the  projectile  to  traverse  its 
diameter,  or 

T0  = 

The  beam  mode  shape  meets  all  special  requirements  of  the  differential 
equation,  independent  of  time.  Therefore,  the  assumed  particular 
solution  is: 

Vx’  fc)  ”  N  -  iVN(x)TN(t) 
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MV 
P  P 


F  *  T  5(x-L)[U(t)-U(t-T0)] 
F  »  ?  6  (x-l)  sin  (irt/T0)[U(t)-U(t-T0)] 


or 


Compatibility:  MpVp  *  /Jp^yU.ydx 


y(0,t)  *  0 

y1  (o,t)  =  o 


A)  Forced  Vibration  At  A  Point 


MV 

6: 


F  -  F{U(x-l+D/2)-U(x-l-D/2)}[U(t)-U(t-Tft)] 


or 


F-  F{U(x-l*D/2)-U(x-L-D/2)}s1nUt/TJ[U(t)-U{t-Tj] 


Compatibility:  MpVp  *  /^pBABy(x,T0)dx 


.  .  y(o.t)  »  o 

B)  Forced  Vibration  Over  An  Area 


Figure  12.  Schematic  of  Forced  Vibration  Models 
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Substituting  this  into  the  beam  equation,  we  have 

N  I  1  <CWX>TN(t>  "  PbVEIN  I  lWx)VC) 


F«(x  -  L)[U(t)  -  U(t  -  T0)]/Eli 


Applying  orthogonality,  as  before,  this  expression  reduces  to, 

6MVTM(t)  "  (pbV/ei)tm(,:)  c  <*VL>/El*>IU<t> 


u(t  -  tq)3 


For  all  time  except  between  t  =  0  and  t  =  Tq,  the  right  hand  side  of 
Equation  18  is  zero,  thus  T^t)  is  zero  for  the  same  time.  For 
0  <  t  <  Tfl,  the  right  hand  side  of  Equation  18  is  a  constant,  therefore, 
TM( t)  is  a  constant  and  the  equation  reduces  to 

e/Si  c  *VL)/EU 

which  further  reduces  to  an  expression  for  CM, 

VV)/pbW  <19> 

This  is  the  particular  solution  to  the  differential  equation  which,  when 
combined  with  the  homogeneous  solution,  leads  to  the  total  solution  to 
the  forced  vibration  problem, 

yT(x,  t)  -  N  |  jtAj^sinCujjC)  +  B^osCi^t)  +  CN1 4>N(x)  (20) 

Evaluating  the  initial  conditions,  y(x,  0)  s  y(x,  0)  =  0,  leads  to 


yT(x,  t)  -  N  Z  ^[1  -  cos(aiNt)]<}iN(x)  (21) 

This  solution  Is  valid  for  0  <_  t  <  TQ,  at  which  time  the  force  is  removed 
and  the  problem  becomes  free  vibration  with  deflection  and  velocity  Initial 
conditions  from  the  forced  response.  Evaluating  and  applying  these  Initial 
conditions  results  in  two  equations  and  two  unknowns  which  provides  ex¬ 
pressions  for  A„  and  Bu.  The  solution  for  t  >  Trt  is  then. 
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yF<x»  O  =  S  jCj^{sin(u^TQ) sin(iOjjt)  -  (1  -  cos(WjjTq)] 
cosCu^OH^x)  (22) 


The  only  unknown  at  this  point  is  F,  which  can  be  determined  from  the 
conservation  of  linear  momentum  assumption.  Projectile  and  beam 
momentum  are  formulated  in  the  same  manner  as  in  the  Initial  velocity 
analysis,  with  the  additional  requirement  that  the  transfer  is  complete 
at  time  Tq.  From  this,  f  is  determined  to  be 

F  *=  M  V  1/ cl 
P  P  S 

where 

CS  “  N  I  i(2*N<L)KN/<0N6N^Sln((VV 


Hence,  is 

s  ■  VpV1,/v»cs 

The  forced  vibration  solution  becomes 

yT(x,  t)  *  (MpVpMgCs^  £  ^(D/wJtl  -  cos(a>Nt))4>N(x) 

(24) 

for  0  <  t  i  Tq 

and  the  free  vibration  solution  becomes 

*F(x’  t)  -  (MpVp/MjC^  1  1(*H(L)/«J){sin(»HT0)8in(»Ht)  - 

(25) 

[1  -  cos(w^Tg)  ]cos(u^t))4,jj(x)  for  t  1  TQ 

This  solution  can  also  be  modified,  as  was  the  Initial  velocity 
solution,  to  spread  the  force  over  an  area  as  wide  as  the  projectile, 
as  shown  In  Figure  12B.  The  analysis  procedure  is  the  same  as  before 
and  results  in  the  same  solutions  as  Equations  24  and  25  with  two 
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substitutions.  ^(L)  is  replaced  by  4^(1.)  =  []  *  qP/ \  and  is 

replaced  by  =  N  f  1  (2^(L)KN/iDN6Ni)sinuNT0> 


A  second  set  of  solutions  can  be  obtained  by  assuming  the  forcing 
function  to  be  a  half-sine  wave,  as  opposed  to  a  square  wave,  acting 
over  the  impact  duration.  This  assumption  makes  the  force  take  the  form, 
also  shown  schematically  in  Figure  12A, 

F(x,  t)  =  F«(x  -  L)sin(*t/T0)[U(t)  -  U(t  -  TQ) ] 
and  Equation  18  becomes 

-  (pBV/E1)TM(t)  *  [F*M(L)/Eli]sin(Trt/T0)  (26) 

Letting  TM(t)  =  R^sinfirt/Tg) ,  differentiating,  and  substituting  into 
Equation  26  results  in  the  equation 

EI6X  +  <*2pbV/to)rm  *  ?V1)/£ 

which,  when  reduced  and  solved  for  R^,  yields 

*M  *  E*M(L)/MB(“m  +  1’2/to)£  (27) 

This  expression  is  substituted  In  Equation  20  and  solved  with  initial 
conditions  for  AN  and  to  give  the  total  forced  response  as 

yF(x,  O  “  N  I  ^(sinfwt/Tg)  -  (Tt/T0uN)sin(uNt)]^N(x)  (28) 


As  before.  Equation  28  is  used  to  develop  initial  conditions  at  t  =  Tq 
for  the  free  vibration  solution  for  t  1  Tq.  This  solution  Is 

yF(x>  C)  "  N  1  1(-*VTOuN)[1  +  cosUNT0)sin(o,Nt)  - 


sin(uNT0)cos(uNt)]<iiN(x) 


(29) 
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The  conservation  of  momentum  yields 

F  =  (M  D  U irCc) 

P  P  S 

where 

S  -  N  I  lV.111"  + 

Thus,  Rn  is  given  as 

VpV^W^  +  ^O*  (3Q) 

which,  in  turn,  results  in  the  total  response  for  0  £  t  <  Tq  being 
yp(x,  t)  =  (MpDp/1rMBCs)N  |  ^^(LXo.2  +  ir2lj)  ][ain(irt/T0)  - 

(31) 

(’r/T0“N)s±n(u)Nt)  1  VX) 
and  the  response  for  t  >  Tq  being 

y^x,  t)  «  |  xHNa >/v4  +  1,2/T0)1 

{[1  +  cos(ojnT0)  ]sin(u>Nt)  -  sinfu^Tg) cos (a^t)  J^ix) 

A  solution  for  the  half-sine  wave  spread  over  an  area  as  depicted  in 
Figure  128  can  be  obtained  by  making  the  substitution  of 
L  +  Dd/2 

/L  _  <r°r  '*n  the  aPPr°Prlate  Places  in  Equations  31  and  32. 

All  the  aforementioned  solutions  for  forced  vibration  at  a  point 
(Equations  24,  25,  31,  32)  and  also  those  for  forced  vibration  over  an 
area  were  converted  to  computer  programs  for  evaluation  against  test 
data. 
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SECTION  IV 

NUMERICAL  ANALYSIS  AND  RESULTS 

Each  theoretical  solution  derived  was  programmed  for  analysis  on 
the  digital  computer.  The  program  was  written  for  use  at  an  interactive 
graphics  terminal  so  that  several  parameters  in  any  particular  problem 
could  be  easily  varied  and  results  compared.  The  required  input  included 
beam  parameters  (dimensions,  material  properties),  projectile  parameters 
(diameter,  density,  velocity),  the  number  of  modes  to  be  used  in  the 
summation,  and  the  time  period  over  which  to  make  strain  calculations. 
Damping  effects  could  be  included,  at  the  user's  option.  When  included, 
damping  was  modeled  by  decaying  exponential  functions  applied  to  each 
mode.  The  values  of  £  Included  in  the  program  were  based  on  actual 
damping  measurements  taken  from  each  test  specimen. 

Initially,  numerical  difficulty  was  encountered  evaluating  sinh  and 
cosh  functions  for  large  values  of  g.  The  problem  arose  in  calculating 
K  in  Equation  3a  and  stemmed  from  taking  small  differences  of  large 
numbers.  As  g  increased,  so  did  the  values  of  sinh  and  cosh  until  the 
accuracy  of  the  computer  was  exceeded.  The  problem  was  solved  by 
expressing  hyperbolic  functions  in  the  mode  shape  equation  as  exponentials. 
A  complete  derivation  of  the  expressions  used  to  calculate  <t>(x)  is  included 
in  Appendix  D.  Once  the  exponential  equations  were  included,  all  numerical 
problems  were  corrected  and  accurate  solutions  were  obtained  from  the 
program,  which  are  listed  in  Appendix  E. 

Results  calculated  and  included  in  the  output  were  the  maximum 
absolute  value  of  strain  calculated  over  the  specified  time  period,  the 
time  at  which  it  occurred,  the  resonant  frequency  of  each  mode  used  in 
the  summation,  and,  if  desired,  a  plot  of  strain  versus  time.  Figure  13 
is  an  example  of  the  program  input  and  output  and  Figure  14  exemplifies 
the  type  of  strain/time  plot  generated.  The  program  was  used  to  calculate 
theoretical  impact  responses  for  comparison  with  each  other  and  with 
experimental  results. 
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Table  2  contains  values  of  maximum  strain  calculated  by  each  model 
at  the  root  of  a  particular  beam.  The  number  of  modes  used  in  each 
model  was  varied  to  determine  how  many  were  necessary  for  the  strain 
value  to  stabilize.  It  is  apparent  that  the  point  Initial  velocity 
solution  does  not  stabilize,  i.e.,  the  more  modes  included  in  the 
solution,  the  higher  the  calculated  value  of  strain.  Spreading  the 
initial  velocity  over  an  area  decreased  the  rate  at  which  strain 
increased,  but  did  not  stabilize  the  solution.  Both  forced  vibration 
models  which  treat  the  force  as  a  square  wave  stabilize  at  five  modes, 
with  addition  of  more  changing  the  max  strain  value  less  than  4%.  The 
half-sine  wave  models  also  stabilize  at  five  modes  with  less  than  2% 
variation  in  the  max  strain  value  thereafter.  Table  3  shows  the  effect 
of  damping  on  the  maximum  strain  values.  The  values  calculated  by  the 
point  and  area  initial  condition  models  decrease  by  24%  and  19%, 
respectively.  Values  of  strain  from  both  square  wave  forced  vibration 
models  decreased  2%  while  the  half-sine  wave  models  predicted  strain 
6%  lower  than  those  with  no  damping. 

These  results  lead  to  two  conclusions  concerning  the  theoretical 
models.  First,  the  initial  velocity  solutions,  although  shown  in 
Reference  3  to  simulate  hard-body  impacts,  do  not  provide  reasonable 
simulations  of  the  soft-body  problem.  The  non-converging  strain  values 
calculated  by  both  are  not  realistic.  For  this  reason,  these  models 
will  not  be  considered  in  test  data  comparisons.  Second,  all  forced 
vibration  solutions  appear  to  provide  reasonable  results  and  should 
be  further  evaluated  and  compared  to  test  data  to  verify  the  validity 
of  each. 


f 
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TABLE  2 

THEORETICAL  RESULTS  WITHOUT  DAMPING* 

5.6  INCH  BEAM;  pp  =  0.032  1bm/in3;  Dp  =  0.5;  Vp  =  340  f/s 


PIC  =  Initial  Velocity  at  a  Point 
AIC  =  Initial  Velocity  Over  an  Area 
FV  «*  Forced  Vibration  at  a  Point;  Step-Function  Model 
FA  «*  Forced  Vibration  Over  an  Area;  Stop-Function  Model 
FV'  *  Forced  Vibration  at  a  Point;  Half-Sine  Wave  Model 
FA'  *  Forced  Vibration  Over  an  Area;  Half-Sine  Wave  Model 


*  Strains  are  absolute  value 
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SECTION  V 

COMPARISON  OF  EXPERIMENTAL  AND  THEORETICAL  RESULTS 

Table  4  shows  the  maximum  root  strain  from  one  of  each  size  test 
beam  and  the  predicted  values  from  each  forced  vibration  model.  The 
percentage  differences  between  the  predicted  and  experimental  values 
indicate  that  forced  vibration  with  a  half-sine  wave  force  over  an 
area  provides  the  best  correlation  without  damping  effects  included, 
although  the  half-sine  wave  force  at  a  point  is  only  slightly  different. 
Predicted  results  for  the  5.6  inch  beam  are  significantly  higher  than 
the  experimental  ones  due  to  damping  effects  in  the  third  and  fourth 
modes  which  were  not  accounted  for  in  the  analytical  predictions.  If 
some  reasonable  values  of  damping  are  included  in  the  predictions,  the 
percent  error  is  decreased  to  that  shown  in  Table  5. 

The  results  in  Table  2  indicate  that  in  all  the  forced  vibration 
models,  the  entire  response  can  be  simulated  using  only  the  first  four 
modes.  This  agrees  with  the  Fourier  transform  plots  in  Figure  9. 
Predictions  of  strain  from  the  same  models  show  that  this  is  also  true 
at  75%  span,  which  agrees  with  Figure  10. 

The  differences  between  predictions  applying  the  force  at  a  point 
and  over  an  area  are  less  than  2%.  This  difference  is  small  enough 
that  results  from  either  model  are  acceptable. 

Figures  15A  and  15B  are  plots  of  predicted  strain  response  from 
both  the  step  function  and  half-sine  wave  models  at  the  5.6  inch  beam 
root.  These  can  be  compared  to  Figure  4A.  Both  models  predict  the 
same  basic  response  shape  as  that  determined  from  experiments  over  the 
50  millisecond  interval. 

Comparison  of  Figures  16A  and  17A  to  the  experimental  short  time 
response  in  Figure  6A  shows  that  the  step  function  model  predicts  beam 
response  more  accurately  than  the  half  sine  wave  model  does.  Further 
comparison  of  Figures  16B  through  16D  to  Figures  6B  through  6D  recon¬ 
firms  the  step  function  model  accuracy.  For  final  verification,  Figures 
18A  through  18D  can  be  compared  to  Figures  7A  through  7D,  which  firmly 
verifies  this  model. 
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Investigation  of  the  half-sine  wave  model  shows  that  It  predicts  the 
same  basic  response  shape  as  the  experiments,  but  the  higher  frequency 
(3rd  and  4th)  modes  are  apparently  not  accounted  for  correctly.  Response 
plots  from  this  model  are  shown  In  Figures  17A  through  17D  for  comparison 
with  Figures  6A  through  6D. 


TABLE  4 

COMPARISON  OF  EXPERIMENTAL  AND  THEORETICAL  RESULTS* 


Beam  Length, 
in. ; 

Pp,lbm/in^; 

Vp,f/s 

Experi¬ 

mental 

FV 

FA 

FV' 

FA' 

£root’* 

er,Z 

a 

er,* 

a 

er,Z 

_ 

a 

er,Z 

a 

5.6; .031;370 

0.278 

0.368 

32.4 

0.360 

17.3 

0.322 

15.8 

8.4; .0307;344 

0.231 

0.272 

15.6 

0.250 

8.2 

0.247 

6.9 

11.2; .0303;366 

0.282 

0.353 

25.2 

1 

0.346 

.. 

22.7 

0.312 

10.6 

0.308 

9.2 

14.0; e033;317 

0.254 

0.302 

18.9 

_ 

0.297 

16.9 

0.272 

7.1 

0.269 

5.9 

*  Strains  are  absolute  value;  theoretical  results  do  not  include  damping 


TABLE  5 

COMPARISON  OF  EXPERIMENTAL  AND  THEORETICAL  RESULTS  WITH  DAMPING 


Beam  Length,  in; 
pp,lbm/ in  ; 
Vp,f/s 

Experi¬ 

mental 

FV 

FA 

FV 

t 

FA’ 

eroot** 

a 

WM 

wm 

2QH 

D 

a 

m 

crf* 

a 

5.6;0.031;370 

0.278 

0.303 

m 

0.297 

6.8 

0.252 

a 

0.250 

-10.0 

VOLTS 
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e,  VOLTS 


TIME,  MS 

A}  Root  Location 


TIME,  MS 
B)  25%  Span 


Figure  16.  Strain  In  14.0  Inch  Beam  From  Step  Function 
Model  (5.0  ms.) 
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TIME,  MS 
C)  50%  Span 


VOLTS  e,  VOLTS 


TIME,  MS 

A)  Root  Location 


TIME,  MS 
B)  25%  Span 


Figure  17.  Strain  in  14.0  Inch  Beam  From  Half-Sine  Wave 
Model  (5.0  ms. ) 
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TIME,  MS 

A)  Root  Location 


10'  VOLTS  e,  10'1  VOLTS 
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TIME,  MS 
C)  50%  SPAN 


D)  75%  SPAN 
Figure  18.  (Cont'd) 
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SECTION  VI 

DISCUSSION  AND  CONCLUSION 

The  experimental  data  provided  an  excellent  base  from  which  to 
evaluate  theoretical  models.  Both  the  FFT  system  and  the  high  speed 
transient  recorder  proved  to  be  very  valuable  tools  for  studying  the 
impact  problem.  Fourier  transforms  of  test  data  made  it  immediately 
apparent  that  the  initial  velocity  solutions  were  not  adequate  models. 
Also,  the  very  short-time  response  plots  from  the  transient  recorder 
dispelled  any  questions  of  high  frequency  bending  waves  causing 
responses  the  FFT  system  could  not  record. 

Both  the  strain  values  and  graphic  displays  of  beam  response  data 
allowed  complete  verification  of  the  step  function  forced  vibration 
model.  Correlation  between  actual  data  plots  and  the  predicted 
response  was  surprising  and  proved  to  be  the  discerning  factor 
between  the  step  function  and  half-sine  wave  formulations.  The  only 
question  that  might  arise  would  be  why  the  half-sine  wave  model  more 
closely  predicts  maximum  strain.  The  answer  to  this  may  be  in  the  damp¬ 
ing  present  in  the  beam/fixture  system.  Accurate  modal  damping  values 
could  not  be  obtained  while  the  beams  were  mounted  on  the  impact  range; 
thus,  any  values  employed  in  the  models  were  approximate.  Conceivably, 
if  the  actual  £  values  for  each  mode  could  be  modeled  into  the  computer 
simulation,  values  of  strain  predicted  would  be  much  closer  to  experi¬ 
mental  values.  Additional  damping  could  also  enter  the  system  through 
the  projectile/beam  Interface.  Plastic  flow  and  many  other  undefined 
phenomena  taking  place  in  the  projectile  during  impact  could  cause 
small  reductions  in  beam  response. 

The  beams  tested  were  specifically  chosen  to  approximate  the 
stiffness  of  typical  jet  engine  fan  blades.  Timoshenko  effects  were 
negligible  in  this  base  because  high  order  modes  did  not  lend  any 
significant  contributions  to  beam  response.  For  beams  of  different 
geometric  properties  where  A/t  (modal  wave  length/beam  thickness)  reaches 
ten  in  the  first  several  modes,  the  Impact  problem  requires  Timoshenko 
beam  theory  to  account  for  shear  and  rotary  inertia.  However,  as  long 
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as  X/t  is  greater  than  ten,  an  Euler-Bernoulli  beam  theory  formulation, 
which  models  soft-body  Impact  as  a  forced  vibration,  treats  the  impact 
as  a  step  function  force  In  time  and  assumes  a  complete  momentum  transfer 
from  projectile  to  beam  will  provide  a  good  solution  for  elastic  beam 
response.  The  accuracy  of  results  from  such  a  model  is  good  for  short- 
time  response  and  can  be  improved  for  long-time  response  by  addition  of 
damping  effects.  Future  work  in  this  area  can  build  upon  this  basis  and 
the  fact  that  linearly  scaled  beams  and  projectiles  will  produce  equal 
strains.  This  means  that  the  response  of  very  large  beam-like  structures 
to  soft-body  impacts  can  be  determined  from  analysis  or  testing  of 
scaled-down  models  and  the  results  will  not  require  scaling. 
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APPENDIX  A  -  ESTIMATION  OF  PROJECTILE  VELOCITY 

According  to  Reference  6,  the  stress  at  the  beam  root  is 

o  -  KAXh /21- 
root  0 

where 

X  =  Impact  Location 

h  =  Beam  Thickness 

Ip  =  Minimum  Area  Moment  of  Inertia 


The  values  of  X,  h,  and  Iq  are  determined  directly  from  the  particular 
beam  geometry.  K  is  the  equivalent  beam  stiffness  determined  from 
Figure  5,  Reference  6  as  a  function  of  impact  site.  A  is  the  effective 
first  mode  amplitude  at  the  impact  site.  For  a  strain  of  0.25%,  °root 
is  25,000  psi  which  gives  A  as 

A  ■=  (50,000) IQ/XKh 

For  the  8.4  inch  beam, 

A  ■  (50,000)  (0.0017086)/ (6. 3)  (202)  (.225)  =  C.299  in. 
Assuming  simusoidal  motion, 

V  =  Au> 


where  u>  is  the  effective  first  mode  frequency,  determined  from 

w  =  K/M 


M  is  obtained  from  Figure  6,  Reference  6.  Thus,  we  have 
w  -  202  x  386/. 227  -  586  HZ 


Thus,  the  velocity  Is 

V  ■  Au  “  175  in/s 

Assuming  complete  momentum  transfer  from  the  projectile  to  the  beam, 
the  projectile  velocity  can  be  approximated  by 


v 

p 


VBM/Mp 
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3 

The  projectile  material  density  Is  typically  0.033  lbm/ln  ,  and  the 
diameter  for  this  case  is  0.7  Inches,  therefore, 

Vp  =  (175)  x  (0. 227) /(0. 00584)  -  6800  in/s 

or 

V  =  570  f/s 
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APPENDIX  B  -  STRAIN  IN  LINEARLY  SCALED  BEAMS 


The  elementary  formula  for  bending  strain  in  a  beam  is: 

e  ■>  Mc/IE  M  =  Bending  Moment 

c  -  Distance  from  Neutral 
Bending  Axis 

I  =  Moment  of  Inertia  of 
Beam  Cross-Section 

E  E  Young's  Modulus 


For  two  linearly  scaled  beams, 

I  =  wt3/12 


I2  =  (K)w(K3)t3/12  =  K Ul1 


and 

Therefore, 


c2  '  KC1 


*1  "  "iV1!* 


and  a  ■, 

e2  =  M2Kc  /K^IjE  =  (1/Kj)M2c1/I-1E 


The  bending  moment,  M,  is  created  by  a  force  that  is  proportional  to  the 
momentum  of  the  projectile. 

T 

L  ®Fdt  =  m  V  "*■  F  e  m  V  /T_  m  =  Projectile  Mass 

0  pp  ppO  p  J 

Vp  =  Projectile  Velocity 
Tq  E  Duration  of  Impact 
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For  linearly  scaled  projectiles  and  equal  projectile  velocities, 


3 

m  ,  ■  pirD,/8 
pi  1 

and 

a  0  *  K3pitD3/8  =  K3m  , 
p2  1  pi 

\ "  "A 

and 

T02  "  K(D!/Vp)  “  KT0l 

Therefore, 

F,  =  m  ,V  /T. 

1  pip  Oj 

and 

F2  *  ■‘Sl’pAj  '  k2f 

The  moment  arm  in  each  case,  L, 

is  also 

linearly  scaled  so  that  the 

bending  moments  are 

w  L1 

and 

M2  =  K  FjKL^  = 

The  bending  strain  is  now 

el  *  Mlcl/:[1E 

and 

e2  =  (1/K3) (K3M1c1/I1E) 

This  shows  that  if  the  force  created  is  proportional  to  the  projectile 
momentum,  equal  strains  will  be  created  in  linearly  scaled  beams  impacted 
by  linearly  scaled  projectiles  (which  Is  in  agreement  with  experimental 
results  observed;. 
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APPENDIX  C  -  EVALUATION  OF  TIMOSHENKO  BEAM  THEORY  EFFECTS 

When  modes  of  vibration  are  considered  in  an  analysis,  the  Euler- 
Bernoulli  Beam  Theory  provides  adequate  eigenvalue  predictions  up  to 
a  certain  point.  When  the  effective  wavelength  of  a  particular  mode 
is  of  the  order  of  ten  times  the  thickness  of  the  beam  being  considered, 
shear  effects  start  to  become  significant  and  the  Timoshenko  Theory  is 
necessary  to  make  adequate  eigenvalue  calculations.  Figure  C-l  is  a  plot 
of  effective  wavelength  to  thickness  ratio  versus  mode  number  for  a 
cantilever  beam.  This  Indicates  that  for  the  beams  used  in  the  impact 
experiments,  shear  effects  will  be  negligible  until  the  twentieth  mode 
of  the  smallest  beam.  The  lowest  curve  on  Figure  16  represents  a  beam 
of  thickness  equal  to  one-fourth  its  length.  Analysis  of  this  class  of 
geometries  should  employ  the  Timoshenko  Theory  to  predict  impact  response 
accurately.  Since  only  the  first  four  or  five  modes  were  apparent  in 
test  data  and  the  same  number  were  necessary  in  theoretical  predictions, 
the  Euler-Bernoull i  Theory  is  adequate  to  model  the  soft-body  impact 
problem  as  considered  in  this  analysis. 
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APPENDIX  D  -  NUMERICAL  STABILITY 
♦N(x)  -  cosh(BNx)  -  cos(B^x)  -  {[sinhO^i)  -  slnCB^f)]/ 
[co8h(8(J'£)  +  cos )  1 J t s inh( BjjX)  -  sin(BNx)] 

To  evaluate  this  expression  on  a  computer,  the  hyperbolic  functions 
must  be  expressed  as  exponentials.  This  leads  to 

♦„(x)  «  (l/ae^d  -  K^)  +  (l/2)e_eNx(l  +  V  -  cos(BNx)  + 
KgSin(B^x) 


where 

Kj,  -  _  e‘BNi/2  -  sin(BN£)  1 /[e^  +  e‘^/2)  +  cos(B^)] 

For  =  17.2787595323  (mode  6),  ^  *  1.000000063.  For  a  16  bit  word, 
the  computer  carries  seven  significant  figures,  and  thus,  would  make 
1  -  Kg  Identically  zero.  e^N  for  mode  6  Is  on  the  order  of  3  x  107  so 
that  (1  -  Kg)e^A  would  have  a  significant  contribution  to  sine  and 
cosine  functions.  Double  precision  will  not  change  with  phenomena  for 
higher  modes,  so  a  convergent  solution  must  be  sought.  First,  in  the 
expression  for  K^,  let  cos(SN*)  =  -  l/cosh(&N^)  and  use  the  identities 
relating  cosh  and  sinh  to  obtain 

K„  «  U  -  (2e’  BNi8in(BN£)/l  -  e"2BN*)Hl  +  e'2^/(l  -  e'2^*)  ] 

Substituting  the  series  expansion 

l/i  -  u  -  l  +  u  +  u2  +  u3  +  .  .  . 


Y  -  2<e~2BN*  +  e"4BN*  +  e-6^  +  .  .  . 
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and  letting 
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we  have 

Kj,  -  (1  +  Y)  -  (  1  +  Y){2e-0N*sin(6/)/l  -  e"2V] 

Now,  defining 

CN  -  1  -  1^  -  [2e'e^sin(e/)/l  -  e~2SN*](l  +  Y)  -  Y 

-2eNi 

and  expanding  1/1  -  e  in  a  series,  we  obtain 

CN  -  2e“BNfsin<BN£)(l  +  3Y/2  +  Y2/2)  -  Y 

This  expression,  along  with  that  for  y,  can  be  substituted  into  the  mode 
shape  expression,  <pN(x) ,  to  yield  this  equation 

♦  (x)  *=  e(0NX  ~  0N^sin(8  £)(1  +  3Y/2  +  Y2/2)  -  (Y/2)e&NX  + 

(1  -  CN/2)e_0NX  -  cos(Bnx)  +  (1  -  CN)sin(BNx) 

This  equation  was  programmed  and  values  of  4>N(x)  were  compared  to  those 
tabulated  in  Reference  10.  The  comparison  showed  that  this  form 
provided  a  stable  solution  for  <|>N(x)  for  the  first  20  modes. 
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APPENDIX  E  -  LISTING  OF  COMPUTER  PROGRAM 


1M» 

lie* 

120* 

»mvi9e> 

130* 

148*C 

1S0*C 

lie* 

no* 

188* 

19C* 

200* 

210* 

220* 

233* 

240* 

2S0* 

2b0* 

270* 

2S0* 

250* 

223* 

210* 

320* 

230* 

340. 

250* 

3«0. 

370* 

3S0* 

330*152 

LOCITv 

4C3- 

EA 

410* 

ION  OUE 
420* 


OUERtAYl IMPACT, 5,0 » 

PROSPAM  IMPDAIlNajT.OuTRUT) 

CC/fflGN  D.KP,V,8L.6J.T.e.K.1,S?<»«l,$PAM.IlM>.g«M!.A!*e».au 

I,S.Tlfl,V(20),EPStl0f8  >,TTIB 

*1*  THIS  PROGRAM  CALCULATES  THE  STRAIN  IN  A  CAHTlliCCft  Ktf 
UNDER  IMPACT.  1TPACT  CAN  El  TREATED  MS  I.C.  Oft  FORCED  VISE. 
BNM*aLL*AL*  AM*  A<*  AA*>»XX.1HN 
B< 1 1*1 .87510407**2. 

8(21*4. 63 403 113**2. 

B(  3 '*7.85475744**2 . 

8<41* 10. 925540744*2. 

B(S)*14.1371£8:>94«2. 

8(61*17. 2737535323**2. 

B  ( 7  1  *20 .  -12035225 1 3**2. 

E(31*23.5£13445C13**2. 

B  < S  '•26.70253755*  ?2. 

El  101*23. 8451333  It. *2. 

Bill  i*32. 5S£7££S5*T2. 

B( 12>*25. 12831551**2. 

B l 13)*39.2iS5C9l7}*2. 

BU41-42. 41150082**2. 

B ( 15  1*45 .5S20934 7**2. 

Bi 161*43. £3463513**2. 

El  17  * *S 1 .23527373**2. 

BUS  1*S  4. 97737 14  3*  *2. 

B(19)*53. 11946403*12. 

B<201*tl. 26105674**2. 

JK*0 

PRINT  152 

FORMAKX. ‘NODEl  OPTIONS  ARE-*, //.10X, ‘POINT  I.C.-INITIAL  WE 
1AT  A  POINT  *,/’lOX,*AREA  I.C.-INITIAL  VELOCITY  OVER  AN  AR 
1,/lOX, ’FORCED  VIBRATION  AT  A  POINT 


•,/ieX, 'FORCED  UIDRAT 


...  ,^IOX,  END  ) 

4j0*  178  PPIH7  1S3 

440-1S3  FORMATi X, ‘OPTION-*  ) 

450*  READ  254, AB 

460*  254  FORMAT <  A2  1 

470*  IF(A5.E0.1HE)C0  TO  176 

480*154  FORMAT (All 

PRINT*. -CO  YOU  UANT  DAMPING?-* 

READ  154, AD 
IF(AA.EQ.1HY1G0  TO  201 

203  PRINT*, *HOU  MANY  MODES  DO  YOU  UANT  TO  USE  9  -  * 

READ*, MM 

IFCBNM.EO.IHYICO  TO  202 

199  PRINT*. ‘ENTER  THE  DIMENSIONS,  DENSITY,  AND  MODUtUS  OF  YOUR 


430* 
SCO* 
SlO* 
520* 
530* 
540* 
550- 
BEAM. • 
560* 
570- 
530* 
590- 
OF  L  - 
600- 


25  PRINT*, •L.U.T.RHO, AND  E  -  * 

READ*.CL,EU,T,R,E 
IF1AL.E0. 1MY1CO  TO  E82 

2«5  PRINT*, 'ENTER  IMPACT  POINT  AND  POINT  FOIk STRAIN  CALCS  IN  ft 
1* 


p;G,  ,3  BEJI  JUALIO  PB4<ai'1*£U1 

J^l***^”*0  - 


57 


6 £*• 
S3*. 
C40* 
6£0*  2 
TILE  - 
ECO* 
670- 
660. 
6S0-206 
ER?«5O0 
700* 
710* 
72©. 
730* 
7*0" 
750'C  1 

760-  2< 

770- 
720. 

7'3C« 

£00*  20 
£10- 
823- 
S33« 


SPAM-SPANt/fO*. 

$P*N2**P4K>.-|00. 

IF'ALl.EJ.I**  ISO  ’C  »2 

2M  print*. -enter  Dianna.  MNtjto.  «mc  * loeiTv  or  t*  p*jcc 

i* 

RE»C*.D.RP.v,' 

IF<  AC.E0.1HVIG0  TO  2*1 

6  PRINT*. -UHAT  TIDE  PERIOD  DO  VOU  WANT  TO  CALCULATE  STRAIN  OV 


i.  so,  s,  or  a.s  ns.»  -  • 

PEnlU.Tin 
TTin-rin 
TIN*  TIM* .001 
IF  (PICK'  .EO.  1  MV  >C0  TO  19 
****  CALCULATE  RESONANT  FREQUENCIES 
£02  Cl-S. 67157*T*SQRT(E/R  >/BLt*2.0 
DO  £  i-i.rm 
u> I )*B( I  )*C1 
6  CONTINUE 

201  IF  >  AB . £0 . 2MPI 'CALL  PIC 
IF(AB.£0.2HAI  'CALL  AIC 
IF(  AB  .  EO.ZHF'J  (CALL  FU 
IF ( A3.E0.2HFA 'CALL  Fa 


fi-io-c  mumimmui  calculate  strain  tttnnmmtiuii 


350- 
£60* 
870- 
££?• 
SO  3- 

so?. 

910- 
920- 
930* 
940- 
950. 
960* 
970. 
930- 
930* 
1000- 
1010" 
1020- 
UtJ  UTI 
1030. 
10*0" 
lose- 
1060- 
1070* 
1080* 
10S0- 
1100- 
1110- 
1 120-C 
1130* 
11*0" 
1150. 


T0.D/(U«12.  1 
13  TI-0. 

S-0. 

IT-0. 

00  12  1-1,1000 
TI«TI*.001*Tin 
EPS* 1 1-0. 

2-0. 

po  i3  j-i.nn 

if<j.eo.i.ahd.l.le.9. iz-e.ioi 

IF ( J .CC.2 -OR. J .E0.3  >2*0.0005 
IF  <  J . Ct .  4  )Z*0 .0015 
IFiAD.EQ. 1HN )Z*0, 

IFiAB.E0.2HFy )C0  TO  75 
IF(AB.E0.2HFA)G0  TO  75 
EPS1-V*  J  i*5IN(U(J)XTI  )*EXP(-Z*U(  J  )*TI ) 

TO  13 

75  EPSl-Yt J ) t ( SI N(u< J )*T0 )*SIN(U( J )*TI )-< 1 .-COS(U< JltTO ) )*COS( 
1 )  )*EXP* -Z*U( J  )*TI ) 

13  EIPSU)^PS?IM!p;rJ>*,‘*C0S,U,J,tTn>*E><P<'2,U‘J,*Tn 
EPS<  I  l-EPSII  UT/I. 

IFlA3S(EPS(II).GT.ABSCS))C0  TO  41 
CO  TO  12 
41  S-EFSU) 

IT- 1 

12  CONTINUE 

PRINT  RESULTS  mnmimmmtti 

Tn-IT 

Tn.Tn*Tin*.o0i 

PRINT  300.S.TH 
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1160*  30*  rORf»T(2X,*T>«  %UC  5T**I«  IS  *.£ XCv*®t*6  OT  MW  « 

sec  0 

1170*  1NPS.*> 

ltso-  print*.*  note  r»€»uiHcv.  hz.  coefficient  »:*, 
IN. 

1190*  1VEL,  F  /S * 

1200*  DO  10  N-l.NN 

lcia-  up*u(n>'C2.*3.14158> 

1220*  A1*A(N> 

1250*  PRINT  101,N,UO,A|,D.V 

12a8*  let  F  OFPATtCX,  I2,7X,E12.6,SX,E12.6,2X,C9.3,2X,C9.3> 

1250"  10  CONTINUE 

12t0*  30  FPl'.TI,  ‘PLOT?  -  « 

1270*  PEhP  15-i,AIJ< 

1220*  IF' AIJC.EO. 1HY>C0  TO  31 

1290-  GO  TO  2S 

1300*  31  CALL  PLOT i EPS  1 

1310-  29  PRINT*,  ‘COULD  VOU  LICE  THE  SANE  PROBLEfl.  DIFFERENT  OPTION? 


1320* 

1330- 

1340* 

1350* 

1350* 

1 3  ’0* 
13SC* 
1390* 
H03* 

CT  OR  S 
H 10* 
HOO* 
HOC* 
t-ne- 
CTILC? 

use* 

h;o* 

H  7  0  • 
14S0* 
CTIIE 
1-iSO* 
1500* 
1S10* 
1520* 
CTILE 
1930* 
1S*0* 

1  998* 

I960* 

1970*176 

1923* 

1950* 


1600* 

NNv,90) 

1610* 

1628* 

1638* 

1640* 


REmD  1S4,AA 
IFtAA.EO. tMVlCO  TO  178 

PRINT*,  •'.'COLD  VOU  LICE  TO  CHANGE  THE  NUNBER  OF  N0DES7  -  • 
READ  lSJ.PNH 
IFiEliN.EO.  1HV )G0  TO  203 

PRINT ( .  "UC'ULD  VOU  LIKE  TO  CHANCE  THE  TINE  PERIOD?  -• 

FEAC  1;4,A1,< 

if*akk.eo.ihv;co  TO  206 

27  PRINT*, 'UOULD  VOU  LIKE  THE  SANE  PROBLEN  WITH  DIFFERENT  INPA 


1TPAI1I  LOCATIONS?  -• 

REAP  15J.ALL 

IP > ALL.EO . 1HY ICO  TO  205 

22  PRINT*, ’JOULO  VOU  LIKE  THE  SANE  BEAN  UITH  A  DIFFERENT  PROJE 


1-* 

READ  1S4.A< 

IF(A<.E0.1HV)C0  to  200 

23  FR1NT*, ‘UOULD  VOU  LIKE  A  DIFFERENT  BEAN  UITH  THE  SANE  PROJE 


1?  -* 

READ  154, AL 

IFiAU.EQ. 1HV1C0  TO  199 

24  PRINT*, ’UOULD  VOU  LIKE  A  DIFFERENT  BEAN  AND  DIFFERENT  PROJE 


1? 

READ  154, AN 

IF (AN.E3. IHVJCO  TO  199 

GO  TO  178 

CONTINUE 

END 

SUBROUTINE  PLOT! X  ) 

CONnON  0,  RP,U, 8L, BU.T.R.NN.SPANI ,$PAN2, B< 20  > ,U< 20 ) , A(20 ) ,DU 


i.s.Tin,Y(20),EPS(ieeo),TTin 
DIMENSION  X( 1080), IT<8) 

CALL  IN  ITT (120) 
lFlA8SlS).LT.e.0ei42lS)VN*.5 


^IHS  1'Avi L  lu  :  .  I  v v  \LlIt  1  K.W 
BVC'il  uXU  U  i  .  . .  .  .  A 
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I6?e 

1660 

1679 

1GS0 

1630 

1760 

1710 

1720 

1730 

1740 

1 7S0 
1760 
1776 
17£0 
1750 
1  SCO 
1310 
1330 
1330 
13-10 
1360 
I860 
1870 
1333 
1S90 
1900 
1910 
1920 
1933 
1943 
1950 
1960 
1370 
1930 
1953 
2030 
2010 
2020 
2030 
2040 
2050 
2060 
2070 

2OS0 

c050 

21C0 

2110 

2130 

2130 

2140 

2150 

2160 

2170 

21S0 

21S0 

2230 

2210 


IFCABStSl.CT. 0.00 1*215  »VH*1  • 

If  <AES< Sl.OT. 0.0023*3  IV1*2. 
CALL  CU!ND0C0.,T2H.-VF,VM 

/-all  T3MI  783.2S881 


TI'Tir.I.COl 
CALL  "OUEAI 0 . ,0. 1 
DO  10  1-1,1000 
X 1 -x  C I  1*2. /0. OOS636 
CALL  DRAUACTt.Xl) 

10  ti*t:*7ipu.ooi 
CALL  ROLEACO. ,-VW> 

CALL  DRAUAC  TIP!,  -VII) 

CALL  LRAUAaiPl.VPl) 

CALL  DRAUACO.  ,VPI> 

CALL  DRAUAiO .  ,-VPI ) 

ChLL  TSEND 
R-T1PI/5. 

DO  11  1*1.4 

CALL  rOU£A(Itft,-V11» 

CALL  PRA'JAl  I*P,  -Vn+VPl*.0*l 
CALL  nCL'EAClJR.YP!) 

11  CALL  DP  ALIA  ( I  *R ,  VPl-YfU .  0* ) 
CALL  TsEND 


R*-YP! 


DO  12  1*1,9 

R.P-YPU.2 

Cw'-L  PtOLEAl6..R) 

CALL  DR«L!A(  .021TIf1,R) 

CALL  NOUEAi  TIPI.R  1 
12  CALL  D°i-UAiTin-.02mi1,R) 
CALL  T5EIID 


J-262 

CALL  ANPIOCE 
CALL  ChRSI2( * ) 

DO  13  1*1.6 
CALL  nOL'At?(J,2*0) 

IF ( I  .E'3.1  )IBL*2H0. 
IFCTTIM.EO.SOOIGO  TO  1001 
IFCTTjl1.E0.S0 >G0  TO  1002 
IFCT7in.EO.5)C0  TO  1003 
1F(1.E0.2>LBL*2H.1 
IF i I.E0.31LBL-2H.2 
IFCI. £0.4  HBL-2A.3 
I F < I .  EO  .5  ILBL-2H.* 
IF(I.E9.6)LBL*2H.5 
CO  TO  1100 

1002  1F(2.EO.2)LBL*2H10 
IF ( I .EO. 3 )LBL*2H20 
IFCI .EO. 4 1LBL-2H30 
IFCI.EO.5)LBL*2H*0 
IF ( ! . EO . 6  >LBL*2H50 
CO  TO  1100 

1003  IF ( I . EO . 2  )LBL*2Hl . 
IF(I;E0.3)LBL*2H2. 

IFCI .  F.Q  .4  >LBL*2H3. 
1FCI.EQ.51LBL-2A*. 

IFC I .EQ.6  )LBL*2H5. 


^  t  c;  --  ‘  L  - 

i  r. .  7  -  ’  ■  *- 


la 


v  o'jALH'i  ? 
i'o  vl>C 


BACXlCAW*8' 
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CO  TO  1100 

223#-  1001  IF(!.C0.2)l»t.-£»n# 

22- 10* 

IFd  .CO.  3  <LBL-cH.* 

2250- 

IF i l.£0. 41LBL-8H30 

2££0* 

1F(!.EQ.5.'UI-£H40 

2270- 

IF(I.EO.6>LIL*£HS0 

22S0- 

CALL  AOUTST <  2, LBL • 

2253* 

LEL-2H0. 

2300* 

IFd.EO.UCO  TO  1101 

2310*  1100  CALL  A0UTSTl2,LBL  > 

2320-  1101  J-J-TS0 

2330- 

13  IFiI.EQ.5.AMO.TTin.EO.50e)J* 

2340- 

J-2S0 

2350* 

00  14  1-1,3 

23=0* 

COLL  r:0uo55(220,J) 

2370* 

IF (S.GT. 0.00234 3  ICO  TO  21 

2320* 

IF<S.LT.0.O01421S>CO  TO  22 

2350* 

IFd.EO.l  1LBL-2H-1 

24J0* 

IF ( I .EQ.2 JLBL-HH0. 

£410* 

IF( I  .£0. 3  )LBL*2HI . 

2423- 

<30  TO  25 

2433* 

21  IFd.EO.l  HBL-2H-2 

£440* 

IFi I .EO .2 (LEL-2H0. 

2450* 

IF ( I .EQ .3  /LBL-2H2. 

2450* 

CO  TO  23 

£470* 

22  IFd.EO.l  )LEL-£H-5 

2453* 

IF  C  2 .£0.2 )LBL-2H0. 

2430- 

IF i I .£0.3 ILEL-2H5. 

2S00* 

25  CALL  A0UTSTl2,LBL) 

2510- 

CALL  CHRSIZI3) 

2520- 

14  J-J-I344 

2530- 

CALL  TSEND 

£540- 

CALL  NOl'ABSf  2000, 140) 

£550* 

EO  50  1-1,4 

£553* 

ITU  ) -2HTI 

£570* 

1T(2  )-2MflE 

£530* 

ITI 3)-£H,M 

£590- 

I T  £  4 ) ■ £H5 ■ 

2630*50 

CALL  AOUTSTO,  IT(  I ) ) 

2610* 

CALL  KOUABS(80, 30881 

£626* 

DO  51  1-1,5 

2630* 

ITIl >-2HST 

2640- 

IT  f  3  > -2HRA 

2650* 

IT  1 3 ) *2HIM 

2660- 

IT(4>-2H,U 

£670- 

IT  £  5 )- 1H. 

2680- 

IFI 5. LT. 0.0014215 )IT(5)*2H-1 

£630* 

IF  £  5 . LT ,0 .0014215  >N- 10 

2730-51 

CALL  AOUTST (2, ITI I 1 ) 

2710* 

CALL  MCf.c 

2720* 

CALL  HDCOPY 

2730* 

CALL  CHR5IZC2) 

2740- 

CALL  ERASE 

2750* 

CALL  FINITT 

2760- 

RETURN 

2770- 

END 

2780- 

SUBROUTINE  FV 

-A-  J  .  "  '  '  ^ iv. 
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2796  • 

HKY( SO  < 
2800* 
2S10* 
2S20*C 
2330* 
2SHa- 
2350- 
2Seo>c 

tt 

2370* 

ESSO* 

2330* 

2500- 

2910* 

292e- 

.♦( 1 .-C 
2330- 
23-S0- 
29S0-C 
2330- 
297.1. 
2533* 
£593- 
3003* 
3C13-C 
3033* 
3033* 
3043* 
3053-C 
**** 
£053* 
3073* 
30S3- 
30  33* 
3100* 
l.-C(N) 
3113* 
3120- 
3123* 
3H0* 
3153- 
3160* 
nnv<93> 

3170* 

31?0- 

3130. C 

3200- 

3213* 

32E0* 

3230-C 

St 

32<0- 

3250- 

3250* 

3270- 


COflHCN  &.RP,U.0L,eu.T.».l0I.SWJ«l.S*«»«.l««*MM» 

l.S.Tin,V<23) 

DIMENSION  C120* 

*****  CALCULATE  CONSTANTS 
201  FN.3.M1S9JDM3.IRP/6. 

T0"D.'(U*12.  ) 

EI.8UTTSJ3./12. 

mutmil  CALCULATE  VCNl'S  tSSSSSSStSSSStSSSS  *********** 


CO  9  N-l.Nfl 
B 1 "SORT! B ( N  ) ) 

G-2. » CENPl -2 • >B1 >*EXP<-4.*81  )*£XP(-6. SSJ )*EW»(-0. *D1 » > 
C(l«'.-a.(2.*SIN«0tUEXP(-Bl  )*a.*Sl)/(l.-ExP<-2.t»t»l 


ARC'S ItSPANl 

Vi  N )  -ENPI ARG-B1 llSINIBlltll. 


♦3.»C/2.*Ct«2./2. >-G*EXP(ARC>/2 


lifiva.  >T£XP< -ARC  J-COSl ARC )♦( 1 . _C(N 1  )*SIN< ARC ) 

9  CONTINUE 

tttititutt  calculate  cs  uumnmu 

CP-.O. 

CO  17  I • 1 , NM 

T£«?.Y<I  Hit  .-C<N))tSIN<U<I)*T01/lU(I  >SSQRTUa>>) 

151  FORMAT (£13.6) 

17  CS-CS-TENP 

tiumtiumm  calculate  A'S  mutmimtimmut 

CC  3  1-1, Ml 

At  1  ) •£ .  »PHU'1V(  1  )/(C5*R*BU*BL*T*U(  1  )SS2. ) 

8  CONTINUE 

Ktimttitmm  CALCULATE  COEFFICIENTS  FOR  STRAIN  smtsttf 


DO  11  It*  1, Nil 
Bc-SCPT ( B ( N ) ) 

B 1 -SORT i B  <N !  >*SPAN2 

C -2.il  EXP l -2. *82 )*EXPt-4.*B2)+EXP(-6.*B21+EXP(*8.*B2> ) 
5un.Ex“(Bl-S21*SlN(B2)m.«3.1Cx2..C**2./2.  )-G*EXP(Bl 


1/2.  iSENPt  -B1 ) *005181  )-(I.-C<N))*SIN(Bn 
11  V(N).SUIHA(N)1B[N)/BLSS2. 

RETURN 

END 

SUBROUTINE  PIC 

COMMON  O.RP,U.8L,8U.T.R,NN.SPANl.SPAN2.»(2»).U(B*).A(201.OU 

1,S,TIM,V<20> 

DIMENSION  C ( 20  ) 

*****  CALCULATE  CONSTANTS 
231  PM-3. 14159*0**3. 0IRP/6.0 
BM-BL*BU*T*R 
C2.pn*V*12.0/BM 

***********  CALCULATE  Y(N)'S  ******************  *********** 

DO  9  N-1,NM 
B1 "SORT ( B <N ) ) 

C*2 .*(ExP( -2.*Bl I+EXPC-4.SB1 >+EXP<-6.*Bl )»EXP(-8.*B1 )) 
CtN)*-0.(2.*SlNl8l >*EXP(-B1  )*(1 .+C ) )/( 1 . -EXP<-2. SB1 ) > 
ARC'BUSPANl 
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32ff< 

♦»1.-C 

3300* 

3310* 

33eO*C 

3330- 

3340* 

3350* 

3360* 

3370«C 

33S0* 

3350- 

3-<oe- 

3410-C 

tttt 

3420- 

3430* 

3440- 

3450- 

3430* 

l.-C!N> 

3470. 

j43C* 

3490. 

3500* 

3S10* 

3S2C* 

BAY ( 90 ) 
3530- 
3540- 

3SS3-C 
3S60- 
3S70. 
35S3* 
3553* C 

» 

3630* 

3610* 

3620* 

3530* 

3640* 

3G5e- 

36E0* 

»<-A2)> 
3S70- 
3530* 
3330* C 
37C0* 
3710- 
3720- 
3730- 
3740-C 
37S0. 
3760- 
3770. 
3780-C 

tttt 


V(H)»EXPtAW-»l>«0H*«it  )*£**>♦'*■**  * 

t<N»/2.  )*EXP<-AR0l-C0S<Ae6>*«  t.*C<NU*$IN«A»et 
9  CONTINUE  _  _ _ 

tnmmmtttit  calculate  cs  M*e«****mi*M**M*M* 

05*0. 

CO  4  I"  *  .Nil 

temp.vu  **2.*<i.-c<N>)/SOSTiKin 

4  C5*C5*TEnP  *  _ -  „ 

tttSUKUtmtll  CALCULATE  »‘S 
CO  8  I-l.NPI 

All  ).C2*VlI)/(C5iU(I»> 

*l*J«y*«itj****f  CALCULATE  COEFFICIENTS  FOR  STRAIN  ttttttttt 


DO  11  N-l.HN 
82*S0RT 1 8 (N  > ) 


B1*S03T(BUI)  )*5PAN2 

C.2.*<EXP(-2.*Ba)4E<P<-4.*82l4EXP(-fi.«2»4EXP(-8.*B2l) 
Suri-EXPt Bl-82 )XSIN(B2)t(l .  43.iC/2.*Gt*2./c.  >-CBEXP(8l  >/£.♦( 


1/2. >«EXP(-8l )«C0S(B1 )-Cl .-C(N)  ItSINfBl  l 
11  Y(N ("SuntAIN  )I3<N  1/BLXX2. 

RETURN 

EMC 

SUBROUTINE  AIC 

COMMON  D.PP.U.BL, BU.T,R,Nn,SPANl,SP«N2,B<20)#U(29)»A(20>,DU 


i,s.Tin,v(20) 

DIMENSION  0(20  > 

«*»»  CALCULATE  CONSTANTS 
201  PH-3.141S3tD**3.0iRP/6.e 
BM-bUBUtTlR 
C2«Fn*U«12.0/BN 

mmmu  calculate  v(N)'S  ntmttmmtm  ttxttttxttt 


DO  9  N'l.tlN 
Bl-SORTCBiNJ) 

C-c. * (EXPC-2.XB1 >«EXP< -4.*B1 )*EXP(-€.iBl )*EXP(-8.*81 )  > 
C<n;.-GM2.*SIN(B1  )*EXP<-Bl )*(1.40) )/(l.-EXPI-2.*Bl ) ) 
Al.Bl*(SPANH0.5*D/Bl> 

A2*B 1  * ' SPAN 1-0 .5*0/81 ) 

YIN  )• ( 0 . StC l N  )X < EXP ( A1  )-EXP(A2))-(l.-R.S*C(N>)«(EXP(-Al )-EX 

1-SINtAl  )4S!N(A2)-(l.-C(N))t(C0S(Al  )-COS(A2))>m/Bl 
9  CONTINUE 

calculate  cs  mmtmumimimt 

CS-0. 

DO  4  I-l.Nn 

TEMP-Y C I  )*2. t( 1 .-C(N > >/SQRT(B( 1 ) ) 

4  CS-CS.TEr.P 

imnmmmti  CALCULATE  A'S  mnitltuuntuititti 

DO  8  I>1. NM 

A(I).C2(V(I)/(CS>UU>> 

8  CONTINUE 

immmiumi  calculate  coefficients  for  strain  mmm 
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3790- 
3SC0- 
38 1  A- 
3S20- 
3839* 
l.-C!N« 
3840- 
3850- 
3859- 
3370- 
3S20- 
33SO- 
BHV<S0) 
3300- 
3310- 
3920-C 
3920- 
3340* 
3S50- 
3360'C 
3970- 
3980- 
3930- 
4003- 
4010- 
4020- 
4020- 
P(  -02 ) ) 
4040- 
4050- 
4050-C 
4070- 
4030- 
4090- 
41C0- 
4110. C 
4120- 
4130- 
4140- 
4150-C 
4160- 
4170- 
4180- 
41S0- 
4200- 
l.-CCNJ 
4210- 
4220- 
4230- 
4240- 


90  II  H-t.W 

B2-S0RTI  |(N)  1 

c-itirci*' 

$y.*I-ExP<Bl-B2<tSlN<B2>*CI  .O.0J/2.-C*t2./J.  .-CtixP<B» 

1/2.  )*EXP<-B1  l*CCStBl  )-(l.-CCN>>-Stl««BI » 

It  v(ro.sun<Atrm3(N(/BLtt2. 

RETURN 

END 

SUBROUTINE  FA 

COMMON  D.RP.U.BL.BW.T.R.NH.SPANl, SPANS. *<20>.U<2#».A<20).DU 


1,$,TIH,V<20> 

DIMENSION  0(20 

ttttx  CALCULATE  CONSTANTS  X*«M*M8«m 
201  PM-3.14139*DU3.*RP/6. 

T0-D/1U112.  ) 

Sl.SU*Tt*3./12. 

tttxttxtxxtt  CALCULATE  V(N)'S  »)«»»»» 

DO  3  N-l.NM 
81-SG»T(8<N>> 

C-2. » ( EXP ( -2. *51 )*EXP(-4.*B1 WEXP(-6. *81 )«EXP(-8.f»t » J 

C(N).-G*(2.»SIN<Bt)IEXP<-Bl)*(1.40))/(l.-EXP(-a.*Bl)> 

Al-Blt<SPANl«e.5*0/Bl> 

A2.B1KSPAN1-0.STD/BL) 

Y(N).(0.S*C(N)*(EXP(A1 >-EXP£A2 ) >-( l.-0.S*C<N> 1* (EXP (-At )-EX 

1-SIN(A1 )*SIN<A2>-<1.-C(N)>*CC0S(A1>-C0S<A2>>>*BI/B1 
9  CONTINUE 

itmimmit  calculate  cs  mmmttnmi 

cs-c. 

DO  17  I-l.NN 

TEMP* Yl  I  )*(l.~C(N))tSIN(U(I  )*T0 J/( U( I  )XSORT(B(  I ) ) ) 

17  CS-CS-TEMP 

xxttxttxxxixtxtt  CALCULATE  A'S  xxxxtxxxaxxxxx 
DO  8  I-1.U1 

A( 1 1 -6 . »PfUV*Y( I )/(CS*R*BU*BL8T*U( I >*«. > 

8  CONTINUE 

ttxxxtttxttxt  CALCULATE  COEFFICIENTS  FOR  STRAIN  *««* 

DO  11  N-l.NM 
B2.SQRT(B(N) ) 

B1.S0RT(B(N))*SPAN2 

G"2.*tEXP<-2.*32  l-EXPt -4 .  *B2 >*EXP(-6.XB2 )*EXP(-8.*B2> ) 
SUN-EXP  (Bl-62USIN<B2>m.  *3.  *C/2.*G*»2./2. )-G*EXP(Bl 

1/2. >*EXP(-B1 HCOSCBl )-< l.-C(N)  )*SIN(81 ) 

11  V(N).SunsA<N)*t(N>/8LX*2. 

RETURN 

END 


: 
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OM> 

irvcji 

2238- 

2*<0-e 

2259- 

2209- 

28?9- 

23D0-C 

** 

2859- 
29C9- 
2910- 
2929* 
2932* 
2949- 
•4I1.-C 
2950* 
2989- 
25’9-C 
2920- 
2S90- 
32;o- 
39:9* 
1*2.  ■'TO 
3229* 
3230- 
39«9-C 
3939- 
3909- 
081?.  ) ) 
3270- 
3230-C 
Mil 
3290- 
31C0- 
3110- 
3120- 

3139- 
l.-C(N) 

3140- 
3159- 
3169- 
3179- 


sur*euri*t  fu 

cwjmn 

i.s.Tin.v(2o>.e?s«j*e«>.m»i.* 

31r£<;51CN  C<23> 

XXXI X  CALCULATE  CCKSTHIT* 

291  Ffl-3.l4is:»0l*3.*aP'i. 

T0-D/CU812. ) 

BI-CU1TX83./12. 

lillllcHl  CALCULATE  VINl't  111111111111111111  IIIIIUIII1 

D9  9  N-l.NM 
EI‘S02T(6(N1) 

C-2. 1 lEXPC-2.181  XEXPI -4.181  >*EXP< -6.811  >4EX?(-8.*B1  »» 

Ci;.'  )--04<2.lSIrf(ll  >8EXP<  -II  ):<l.*C)>/(t.-EXP<-2.XXt >> 

arc-?  ixsfam 

V(N>-EXP(ARC-B1>8$1N(B1  >*(1.43.*C/2.4CtT2./2.  >-C8EXP(AftC>/« 

1(N)'2. )*EXP< -ARC )-C0S( ARC  14(1. -C<M) IXSIn(ARS) 

9  CONTINUE 

*8111111118  CALCULATE  CS  X88X888X8XXXt8 

C5-9. 

CO  17  1-1. Nil 
C5-S3RTt8(l ))/3L 

Tcr?.(l.»COS(U(l  )8T0))8(l. -C< I) )8V( I )/(BS8<U<I 1882. «3.14IS9 

1882. 1) 

17  CS-CS4T£nP 

8888111181(181888  CALCULATE  A'S  888*888*888888*888888*8* 

DO  8  I-l.NN 

A( 1  ).P-101V( 1 1/(283. 14 1S98R8SU8T8CS8IUI 1 1882. 43. 141S9882./T 

8  CONTINUE 

*8811X11118888X88  CALCULATE  COEFFICIENTS  FOR  STRAIN  8888X888* 

DO  11  N-l.Nfl 
B2-5GRT(B(N>) 

81-SCRT(B(N))8SPAN2 

C-2.8(EXP(-2.8B2)4EXP(-4.*B2)4CXP(-6.*B2)4EXP(-a.*B2>l 
SUN-EXP (B1-B2>XSIN( B2)*(l -43.tc/2.4CS*2./2. >-6XEXP(Bl 1/2. 4( 

l/2.)8EXP(-B114C0S(9n-<l.-CIN»)*SlN(Bl> 
tl  V(N)-SUntAlNmifi>/BL*82. 

RETURN 

END 


FRACIICiJSLK 


1.130 
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3029- 
my  <901 
MM* 

3340- 
3759 -C 
3050- 
M74- 
3ISB- 
3S00-C 
4032* 
4»ie* 
4920* 
4939* 
4049* 
4050* 
4969* 
P<-A2>> 
49?9* 
40£3* 
4093-C 
4100* 
4110- 
4129* 
4133* 
M2./T0 
4149* 
4153* 
4160-C 
4179- 
4188- 
0*12. ») 
41S0- 
4299-C 

4218* 

4228- 

4239- 

4249- 

4250- 
l.-CCN) 

4269* 

4279- 

4239- 

4290- 


*ub*outjn*  m 

I.I.T19.VIMI 
diiznsScn  ci£ti 

*****  calculate  constants  ***i*t*xat*a* 

291  PN0.14ir;jSt*3.(«P/9. 

T0-D'CV*12. J 
8WUtT3*3./!2. 

************  UlCUKTt  V(Ni*S  ************* 

DO  B  N-I.NN 
*1-SCRT<»(N>> 

C-2.  *C£xPC -2.  *81>»ExPl-4.  til  >*EXPC-9.*81>4tXPC-B.»Bl>> 
CCNI--C*<2.*SINC*t  )*EXPC-»1>*C1.40>>/C1.-£XPC-2.*»I1> 
Al-BI*CSPAN1*9.S*D/BL> 

A2-lUCSPAIil-0.S*0/*L) 

V(N>-C8.S*CCN)*CEXPCA1  J-EXPCA2 ))-(!. -9. S*C<N>>*CEXP<-A1 »-EX 


1 -SIN Cnl  )«SINCA2  >-( 1 .-CCNi XCCOSCAl  1-C0SCA2) ) )*3l/li 
9  CONTINUE 

xtscist******  CALCULATE  CS  ***************** 

CS-8. 

DO  17  1-1. NH 
BS*S02TCBCI))/Bt 

TE“P-  ( 1 .  «C0SCU(  1  »*T0))*tI.-ClI >1*VC  I  )/( BS* CUC I >*X2. *3.14199 


1**2. >> 

IT  CS-CStTECIP 

ztxttxsis******  CALCULATE  A'S  *************** 

00  8  1-1, NX 

AC  I  )-PA*D(V( 1 )/<2*3.141S9tCSsRX3U*T*CUC 1 )**2.»3.141S9**2./T 


8  CONTINUE 

*************  calculate  coefficients  for  STRAIN  ***** 

DO  11  N-l.NH 

B2-S32TC8CN)) 

s: -S03T C BIN )t*SPAN2 

C-2.*CEXPC-2.132MEXPC-4.*82>*EXP«-6.*»2l4£XP«-8.*B2>> 
SUn-EXPCBl-B2!*SlN(B2)*tl.43.*S/2.+G**2./2.1-ttEXPUl  )/*.♦! 


1/2.  )*EXP(-B1 UCOSCBl )-( 1 .-CCN))*SIN(B1 
11  VCN>-SU.-l*ACN>*8CH)/*l**2. 

RETURN 

END 


> 
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